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We consider a scenario where interacting electrons confined in quantum dots (QDs) are either 
too close to be resolved, or we do not wish to apply measurements that resolve them. Then the 
physical observable is an electron spin only (one cannot unambiguously ascribe a spin to a QD) and 
the system state is fully described by the spin-density matrix. Accounting for the spatial degrees of 
freedom, we examine to what extent a Hamiltonian description of the spin-only degrees of freedom 
is valid. We show that as long as there is no coupling between singlet and triplet states this is 
indeed the case, but when there is such a coupling there are open systems effects, i.e., the dynamics 
is non-unitary even without interaction with a true bath. Our primary focus is an investigation of 
non-unitary effects, based on exact master equations we derive for the spin-density matrix in the 
Lindblad and time-convolutionless (TCL) forms, and the implications for quantum computation. In 
particular, we demonstrate that the Heisenberg interaction does not affect the unitary part (apart 
from a Lamb shift) but does affect the non-unitary contributions to time evolution of the spin-density 
matrix. In a sequel paper we present a detailed analysis of an example system of two quantum dots, 
including spin-orbit effects. 

PACS numbers: 



I. INTRODUCTION 



In many quantum computation proposals the spin of 
a localized particle, e.g., an electron or nuclear spin is 
a natural carrier of quantum information - a qubit - 
which can be effectively protected and/or processed to 
achieve a computational task. The ability to govern the 
spin state via controllable interactions is a key ingredient 
underlying several proposed scalable quantum computer 
architectures in semiconductor nanostructures where the 
spin of an electron localized in a quantum dot (QD) or by 
a donor atom serves as a single qubit jiiSi&iiiLSiLSai i 1 1 1 1 It 
is a typical assumption that a single electron is trapped 
in each individual QD and an electron spin can unam- 
biguously be assigned to a QD. Using the spin degree of 
freedom of electrons trapped in QDs (rather than their 
charge) for information processing is of special interest 
since spins have comparatively long coherence times in 
semiconductor nanostructuresj 12 i 13 ' 14 i 15 

In this work we revisit a rather fundamental issue. In 
essence, we ask: "what is a spin?" . To clarify this ques- 
tion (which is not meant in the sense of the relativistic 
origin of spin, as our work is essentially nonrelativistic 
but includes relativistic corrections such as spin-orbit in- 
teraction), we distinguish between the notion of a pure 
spin and a pseudo spin. 

Let us start with pseudo-spins, a notion which applies 
to the majority of studies utilizing spins for quantum in- 
formation processing. In order to ascribe a spin to a local 
site, such as a QD, one defines the electron spin opera- 
tor as a bilinear combination of electron annihilation and 
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(s is a spin index, A is the QD index) 
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(see, e.g., Appendix A in Ref. [l6l ). Then the operators 
{s^ja obey the usual su(2) commutation rules. Eq. ([1]) 
implies implicit dependence of the spin operators on coor- 
dinate degrees of freedom and allows one to ascribe a spin 
to a QD. We call a spin defined by Eq. ([I]) a pseudo-spin, 
meaning that it carries some coordinate dependence. It 
is important to note that the Heisenberg exchange inter- 
action arises only between pseudo-spins. This is a con- 
sequence of the coordinate dependence of pseudo-spins. 
Indeed, the Heisenberg exchange interaction constant Jh 
is electrostatic in nature, i.e., it derives from the Coulomb 
interaction^ 

As opposed to a pseudo-spin, we define a pure-spin as 
a spin that does not have any coordinate dependence. 
This is the usual definition of a spin operator via the 
Pauli matrices a = (cr x ,a y ,<7 z )rl which do not depend 
on coordinates. Considering the spin degrees of freedom 
as carriers of quantum information, the spatial degrees of 
freedom must, in principle, be irrelevant for the storage 
of quantum information. From the quantum information 
point of view it does not matter in what orbital state an 
electron is if one neglects possibly small spin-orbit inter- 
action effects. Although it is not possible, physically, to 
divide the system into purely spin and spatial parts [as 
the "spin" physics is embedded in space], we would like, 
in a succinct description of spin dynamics, to eliminate 
irrelevant information. However, rather than ignoring 
the coordinate dependence, we follow the standard pro- 
cedure of open quantum systems 18 and trace over the 
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coordinates, leaving us with a spin-density matrix as the 
primary object of investigation. Our motivation for in- 
vestigating such a pure-spin model is a scenario where 
interacting and confined electrons either are too close to 
be resolved, or one does not wish to apply measurements 
that resolve them. Then the physical observable is an 
electron spin only (one cannot unambiguously ascribe a 
spin to a site such as a QD) and the system state is fully 
described by the spin-density matrix. 

For simplicity of analysis we consider a prototype 
system of two interacting electrons trapped in two sites 
A and B and separated by the distance tab- One of 
our first findings is that the Heisenberg interaction, 
except for the Lamb energy shift, leaves the unitary 
part of the pure-spin dynamics invariant. This means 
that the pure-spin setting is inappropriate for any of the 
multitude of approaches to quantum computation and 
decoherence control which rely on enacting quantum 
logic gates via control of (sometimes only) Heisenberg 
interaction; 1 - 2 - 3 - 4 - 5 - 6 - 7 - 8 - 9 - 10 - 11 - 19 - 20 - 21 - 22 - 23 - 24 - 25 - 26 - 27 - 28 - 29 - 30 

In these cases one must be able to resolve spins, as in 
the pseudo-spin setting. After developing the appro- 
priate formalism we turn to a comparison of pure and 
pseudo-spin dynamics, in particular spin-orbit effects. 

In the case of pure-spins we are particularly interested 
in finding the conditions under which their dynamics is 
unitary (effectively closed system) or non-unitary (open 
system). We show that as long as there is no coupling 
between singlet and triplet states the dynamics is uni- 
tary (i.e., a Hamiltonian description of the spin-only de- 
grees of freedom is valid), but with the singlet-triplet 
states coupling there are open systems effects, i.e., the 
dynamics is non-unitary already due to the orbital de- 
grees of freedom, even without coupling to a true bath. 
To exhibit these effects we derive several master equa- 
tions for the spin-density matrix, and analyze their im- 
plications for quantum computation. Our central results 
are contained in Eqs. (gSJ) and Eq. is a 

Lindblad-like master equation for the pure-spin dynam- 
ics, which clearly exhibits the non-unitary nature of this 
dynamics (it does not, however, invoke a Markovian ap- 
proximation). Eq. (|M|) does the same within the time- 
convolutionless approach. 

This paper is the first in a series of two. In this pa- 
per (part I), we develop a variety of general models for 
pure-spin open system dynamics. Specifically, in Sec- 
tion HH the operator sum representation for the spin den- 
sity matrix is derived, while the Lindblad-type and time- 
convolutionless (TCL) master equations are considered 
respectively in Sections IIIII and IIVI We conclude Part I 
with a discussion in Section fVl 

In part 11^ we highlight the differences and relation- 
ship between pseudo- and pure-spin models, and provide 
a concrete illustration in terms of a system of two quan- 
tum dots trapping one electron each. In particular, we 
present calculation results demonstrating non-unitary ef- 
fects in the pure-spin model due to both external mag- 
netic field inhomogeneity and spin-orbit interaction. 



Atomic units, H = e = m e = 1, 1/c ~ I/I37, are used 
throughout the paper unless otherwise stated. 

II. THE OPERATOR-SUM REPRESENTATION 

A. Hamiltonian 

The total two-electron system Hamiltonian has the 
generic form 



H = hi + ll2 + Wl2 



(2) 



Here hi, i — 1,2 is a one-electron Pauli Hamiltonian 
which includes spin-dependent terms 

^ = tt-t + -A(n,t)) +V Tr (n,t)+B(f i ,t)-s i (3) 

Zm* \ c J 
where m* is the effective electron mass in the medium, 

31,32,33,34.35.36,37.38,39.40 . . 

{ri,'Pi) &K im eiyi! Irons position and momentum oper- 
ators; A(ri,t) and Vrr(ri,i) are, respectively, the vector 
potential and the trapping potential, which has two min- 
ima at sites ?a and rg where the electrons are localized. 
The magnetic fields 



B(fi,t) = B c ^(ri,t) + B so (fi,Pi,t) 



(4) 



are due to the external (possibly spatially inhomogc- 
neous) magnetic field B cx and the spin-orbit interaction 
field B so which can usually considered to be a small per- 
turbation; Si — \&i is the spin vector of Pauli matrices 
(here we included the gyromagnetic factor g e and Bohr 
magneton hb in the definition of the magnetic fields so 
that magnetic fields are measured in energy units). The 
two-electron interaction term 



W\2 



V cc (ri 2 ) + Vdip(si,S2,ri 2 ) 



(5) 



consists of the interelectron electrostatic interaction po- 
tential V ee {ri 2 ) = ^/{eri2), where e is the dielectric con- 
stant of the medium, with e = 1 in vacuum, and the 
spin-spin magnetic dipole interaction^ 2 - 



y dip = 1.45 meV 



si ■ s 2 r? 2 - 3(s*i • f 12 ){s 2 ■ r 12 ) 



ysi • s 2 S(r 12 



(6) 



which contains a term inversely proportional to the cube 
of interelectron distance r\ 2 , and the contact term pro- 
portional to the (^-function (coordinates and spins in 
Eq. ([6]) are measured in atomic units). 



B. Basis States 

If the system described by the Hamiltonian Eq. © 
can be considered as a closed system that does not in- 
teract with its environment, then dynamics of the state 
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^tot is governed fully by the corresponding Schrodinger 
equation. The state \& t ot {r, o~, t) depends both on the 
electrons' spatial coordinates r = (r 1,7*2) and the spin 
variables a = (cti,(T2). Within the ground state approx- 
imation, which consists of neglecting excited states, we 
assume that an electron can be trapped in two ground 
orbital states 4>A and 4>b localized near the sites ?a and 
tb with energies ea and eb respectively; 4>a and (j> B are 
orthonormal states obtained in the same trapping poten- 
tial. We do not assume a priori any symmetry of the 
trapping potential. With two up and down spin states 
Xm( <7 ) ( s zXl,l — ±3XT,i)i the one-electron basis com- 
prises 4 states: <AA(ri)XT,i( C7 ») and ( pB(n)xt,l( (T i)- Tlic 
corresponding two-electron basis set is defined by Slater 
det erminant s . 1 7 

The singlet subspace: 

*« = fsi(f) X s(o-), i =1,2,3 (7) 
where two-electron symmetric orbitals 

fsi(r) = ^((/>A(ri)(f)B(r2) + (f>A(r2)(t>B(ri)) 

fs2(r) = 4>a(^i)^a(^i), f S 3{r)=<t>B{ri)(j>B{r2) 

represent the states of single and double occupancies re- 
spectively, and the antisymmetric singlet spin function 

Xs = ^ X1 ^ XL ^ ~ ^T( cr 2)xi(o'i)) 

describes the spin state with total spin 5 = and mag- 
netic spin projection number Ms — 0; 
The triplet subspace: 

*« = ft^xtiiv), i = l,2,3 (8) 
where a two-electron antisymmetric orbital 

ft(r) =-= {4>A{ri)4>B{r2) - <t>A{r*2)4>B{ri)) 

and the triplet symmetric spin states 
Xti = Xt( (J i)Xt(o'2), 

Xt2 = -= (xrOikiO^) +Xt(< j 2)xi(c-i)) . 
Xa = Xl(o-i)xi(o- 2 ) 

describe spin states with S = 1 and Ms — 1,0,-1 re- 
spectively. 

Note that the six basis functions $ S j and &u in Eqs. 
([7|) and are orthonormal, as are their orbital and spin 
functions: 

{fsi \fsj) = % (fsi \ft)=0, (9) 

(Xti \xtj) = Sij, (xu \Xs) = 0. 



In the above basis set, we have 

3 

where the expansion coefficients are solutions of the 
Schrodinger equation 

i\a s ) = H ss \a s ) + H st \a t ) , (11) 
t|ot> = H ts \a„) + H u |ot) 

where the column vectors |a s ) = (a s i, a S 2, a S 3) T and 
\a t ) = (a t i, a t 2, at3) T are vectors of singlet and triplet 
amplitudes, and 

Hf? = (<S> si \H\$ sj ), Hg = {* ri \H\* tj ), 

Hj? = (* ti |#|* aJ >, fl*=<$«|#|*y>, (12) 

are correspondingly singlet-singlet, singlet-triplet, 
triplet-singlet, and triplet-triplet subspace interaction 
Hamiltonians. We note that within the Heitler-London 
(HL) approximation^ one neglects the double occu- 
pancy states contribution in the expansion Eq. (|10[) . 

C. Spin Density Matrix 

By analogy with the theory of open quantum systems, 
one can formally consider the spin subsystem as a "sys- 
tem" while the orbital degrees of freedom belong to the 
"bath" . Then the total Hamiltonian Eq. ([2]) governs the 
evolution of this "system + bath", with the "system" 
Hamiltonian being the Zeeman interaction term with a 
space-independent magnetic field (this term may be ab- 
sent), the spin-independent terms describing the bath, 
and the spin-dependent terms (which depend both on 
spin and coordinate variables) describing the interaction 
between system and bath. 

A description of the open system in terms of a com- 
pletely positive (CP) map will result if (but not only 
we assume that system (S) and bath (B) are ini- 
tially decoupled, so that the total initial density matrix 
Ptot{t = 0) = |*tat(0)> (*tot(0)| is a tensor product of 
the system and bath density matrices [p(0) and pb{0) 
respectively], with p(0) being defined on the whole sys- 
tem subspace. The system dynamics is described by the 
spin density matrix: 

p(0) .— » P (t) = Tib [U(t)( P (0) <g> p B (0))rf(t)] . (13) 

Here Tr^ is the partial trace over the coordinates 
( "bath" ) and the time-evolution operator is 

U{t) = TV exp <1tH{t) \ , (14) 

where JV- is a chronological time-ordering operator. 
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Now we come to an important observation: it follows 
from Eqs. J7]) and © that the spin and orbital degrees of 
freedom in |$ tot (0)) [Eq. (fTOf ] are factorized only if the 
state belongs to either the singlet or the triplet subspace, 
but not to a superposition of both. Therefore, the total 
initial density matrix cannot be represented in a tensor 
product form if the initial state contains both singlet and 
triplet parts. In order to understand the case of a CP 
map description, we next consider these two cases. 



2. Triplet initial state 



Similarly, in the triplet case the initial state is specified 



by 



where 



Ptot(O) = p(0) ® |/t) (ft\ 



p(0)=^a ti (0)a| i (0)T ii , 



(19) 



1. Singlet initial state 

In this case |<2t(0)) = and the total density matrix 
state takes the form 

Ptot (0) = p(Q) ® a si(°Kj(°) \f«) (/-J I 



with 



Tij = \xti) (Xtj\ 



(20) 



being coupling operators between triplet states, and 
(at (0)| a t (0)) = 1. Then, the operator sum represen- 
tation is exactly the same as Eq. (|15[) but the Kraus 
operators {A S i,A t } are defined differently as 



where p(0) — S with S = \x s ) (Xs\ being the projection 
operator on the singlet subspace, and (a s (0)| a s (0)) = 1. 
With this initial state, Eq. (fT3|) can be rewritten simi- 
larly to a CP map in the operator sum representation 



as 



18.45.46 



p(t) = A m (t)p(0)Al(t) + Mt)p(0)A\(t), (15) 



where the Kraus operators {A S i,A t } expressed in terms 
of the a s ,t(t) coefficients are 



A sl (t) =a si (t)S, A t (t) = J2a ti (t)Kl (16) 



with 



K > = \Xs) (Xti\ 



(17) 



being operators coupling the triplet and singlet sub- 
spaces. Note that the Kraus operators depend on initial 
conditions via the dependence of \a s (t)) and |a*(t)) on 
the initial amplitudes |a s (0)). This implies that Eq. (JT5J) 
does in fact not represent a CP map - more on this below, 
in subsection III Dl The first term in Eq. (|15p . containing 
the sum over i, and the second one describe, respectively, 
the singlet and triplet states contributions to p(t). 

Using the normalization condition for the a's - 
(a s (t)\ a s (t)} + (at(t)\ at(t)) = 1 - one derives the nor- 
malization condition 



E 

i=l 



AT, A. 



■ A] A t = S, 



(18) 



which guarantees preservation of the trace of p in the 
case of a singlet initial state. 



A si (t) = J2 U ?1 (*) K ^ . A * (*) = E u *j (*) . ( 21 ) 

where the evolution operator matrix elements are defined 
as solutions of the differential matrix equation 



iU = HU, 
17(0) = I 

Ifss JJSt \ / JJSS JJSt 

U = I jjts jjtt J i H = ( jj ts jj tt 



(22) 



Here I is the 6x6 identity matrix. Notice the difference 
between the singlet and triplet cases in the definition of 
the Kraus operators. Namely, the triplet Kraus operators 
do not depend on the initial conditions, unlike the singlet 
case. The initial triplet amplitude dependence is present 
only in p(0). Nevertheless, even in the triplet case we 
do not obtain a CP map, since there is a dependence of 
the Kraus operators on the domain (in this case on the 
subspace of triplet states). 

The normalization condition changes to 



(23) 



where T = Tr Ty = ^ Tu is the projection operator on 
triplet subspace. 



Mixed initial state 



Although a mixed initial state, which contains both 
singlet and triplet parts, cannot be represented in a prod- 
uct form, a slight modification of Kraus operators in 
Eq. (|TB|) according to 



(m) 

si 

M 



(t) = A sl (t)/ (a s (0)\ a s (0)) 1/2 . 
(t) - A t (t)/(a s (0)\ a s (0)) 1/2 



(24) 
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will provide the operator sum representation Eq. (|15[) 
for the mixed initial state. The normalization condition 
Eq. (fT8| is accordingly renormalized so that the right- 
hand side is divided by (<z s (0)| a s (0)). If (a s (0)| a s (0)) = 
1, the mixed case formulas go over into the singlet ones. 



Using Eq. (fT6| for the Kraus operators, the operator 
sum representation Eq. (|15[) can be reduced to 



P (t) = (o.(t) |a,(t)) S+ J2 (M*)> T«. (25) 



D. Analysis 

Note that in spite of having assumed a factorized ini- 
tial state in the singlet and triplet cases, the operator 
sum representation Eq. (|15p . is not completely positive. 18 
This is because the singlet and triplet Kraus operators, 
Eqs. (fT6|) and pTj) . depend upon which initial state (sin- 
glet or triplet) is chosen. This happens since in the 
general case, when the initial spin density matrix p(0) 
contains both singlet and triplet parts, it is not possi- 
ble physically to realize such a state in the product form 
p(0)®\f) (/|, with |/) being a reference coordinate wave- 
function, in the total Hilbert space. In this case, |/) 
would be simultaneously symmetric and antisymmetric 
with respect to permutations of coordinates, which can 
be realized only if / = 0. Therefore, it not possible to 
obtain Kraus operators independent of p(0) in the gen- 
eral mixed case. The Kraus operators in Eq. (|24[) do 
depend on p(0). There is no contradiction to the Kraus 
representation theoreni ) 18 i 45 i 46 which states that a map 
p(0) — > p(t) has the Kraus representation if and only 
if it is linear, completely positive, and trace preserving, 
because there is a subtle difference between the Kraus 
representation of a map and the Kraus representation of 
a state under the action of a map. The Kraus opera- 
tors describing the Kraus representation of a map are 
independent of the state p(0) while the Kraus operators 
describing the Kraus representation of a state p(0) un- 
der the action of a map may be dependent on the initial 
stated 



III. LINDBLAD-TYPE MASTER EQUATION 



While the operator sum representation, Eq. ([15]) . can 
in principle be used to investigate spin dynamics, it is 
not clear how to separate out the unitary evolution of the 
(spin) system from the possibly non-unitary one, which 
results from the system-bath coupling. The reason is that 
in general, each Kraus operator will contain a contribu- 
tion from both the unitary and the non-unitary compo- 
nents of the evolution. In Ref. |48j, a general procedure 
was devised to derive from the operator sum represen- 
tation a CP Markovian master equation, where unitary 
and non-unitary terms could be identified. The Marko- 
vian part of the derivation was based on a coarse-graining 
procedure. Here, we use a similar formal approach, based 
on explicit expressions for the spin density matrix, but 
avoid the coarse graining step since we are not interested 
in the Markovian limit. 



It is straightforward to check that in the general mixed 
case where both triplet and singlet components are 
present in the initial state Eq. (|25|) holds true. After 
some algebraic manipulation one can transform the time 
derivative of Eq. (f2"5| into: 



dp 

at 



[H u ,p] +C ts (o,(t) j0t (t)) 



where 



ij=i 



(26) 



(27) 



and 



C-ts — 



Tr(G)S- VG„T 



ij=i 



(28) 



G, 



F, 



{H ts \a s {t)){at{t)\), 



In the derivation of Eq. (|26|) . we have used Eq. (|lip . 
Clearly, if H ts = (no singlet-triplet coupling), then 
Cts = and Eq. (f2"6"| describes unitary evolution. 

Since the orthogonal projection operators S and T 
commute with p, Eq. (126p is invariant under the trans- 
formation 



H' 



E a (t)S + E t (t)T 



(29) 



where E s {t) and E t (t) are any functions of time. In the 
limit H ts = 0, H" represents an effective spin Hamilto- 
nian. The Hamiltonian Eq. (|29p is clearly Hermitian if 
E s (t) and E t (t) are real functions. 

A procedure to obtain the effective spin Hamiltonian 
was proposed in Ref. [32ll33lf34l . based on a comparison of 
the two expectation values 



(*tot|J?|*tot) =Tr(H spin/ 9). 



(30) 



From the general form of the density matrix Eq. (|25p and 
the relationship Eq. ([30)1 it then follows that H sp i n must 
have the general representation 

3 

H spin = E s (t)S+ E W (*) T y ( 31 ) 

ij=i 

where the functions E s (t), E t ij{t) satisfy the Hermiticity 
conditions: E s {t) = E*(t), E uj (t) = E* tjl {t) . When 
H ts = one derives from Eq. (|3"0|) 



E s {t) 
E ti j (t) 



(a s (t)\H ss \a s (t))/(a s (t) \a s (t)) 



TTtt 



(32) 
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Note that the Hamiltonian Eq. (|29|) fits the general 
representation form Eq. (f3"Tj) . Also, it is clear that the 
Hamiltonian (|3ip is symmetric with respect to permuta- 
tion of spin indices since the basis operators S and are 
symmetric with respect to spin permutations. In terms 
of single-spin operators si and s 2 they take the form 



S 

T22 = 

T33 = 

T12 = 

T13 = 

T = 

T31 = 

where 

Js 

S± 



= ll-si-s 2 , T11 = \l + \S Z + s Xz s 2z 

jl + S lx S 2x + SlyS 2 y ~ Sl z S 2z , 



%T - ±9 

72 12 s - 



SlxS2. 

T\ 3 , 



Sl z S 2 

+ -rJ s ], ±23 = ^L2° + 
c — s ly s 2y + i {si x S 2y + S 2x S\ y 

s*i • s 2 , 



T21 - T 12 , 



-32 



-■t 
-23 



(33) 



si z s 2x + Si x s 2z + i (s lz s 2y + s ly s 2z ) 
S x ± iiSu, 5 = si + S2, 



and I is the 4x4 identity matrix. At the same time, the 
singlet-triplet basis operators Ki are asymmetric 



Ki 



K, 



2^2 



(34) 



where J as = s 2 — Si x S . This symmetry property of 
the Hamiltonian and spin density matrix is quite gen- 
eral since it is a consequence of the orthogonality of the 
coordinate wave functions in the singlet and triplet sub- 
spaces [Eq. ©]. Moreover, it remains valid beyond the 
two-state approximation used here: with inclusion of ex- 
cited states, spatial orbitals lying in singlet and triplet 
subspaces will still be orthogonal to each other. 

We emphasize that the recipe Eq. (|30|) implicitly as- 
sumes that the spin dynamics is unitary. Besides, observe 
that Eq. (|30|) is not invariant under the transformation 
H spin h-> H spin + E s (t)S + E t (t)T. Using Eqs. {33]), one 
can rewrite 



rj/)S - E,(/)T = ( lE,(t) + ^Et(t)jl+ J H {t)s x -s 2 , 

(35) 



where 



J Hit) =Et(t)-E.(t) 



(36) 



and the last term, £/ex, is the familiar Heisenberg ex- 
change interaction. Invariance of Eq. (26]) under the 
transformation Eq. R29\) means that the spin density ma- 
trix will not change under a unitary transformation in- 
duced by this Hamiltonian transformation. In particular, 
it follows that the unitary transformation induced by the 
Heisenberg exchange interaction 



where Aij(i) = L dt' Jn(t') does not affect the state: 

P (0) 1 * p(t) = U H (t) P (0)U H (t) = p(0), (38) 

where p(0) can be taken in the general form p(0) = p s S + 
*52iij PtijTij, with p s and puj being arbitrary parameters 
specifying p(0). It is easy to verify Eq. (f38|) directly using 
the identities S 2 = S, T 2 = T, ST = TS = 0. 

The conclusion that the Heisenberg interaction does 
not affect pure-spin dynamics has important implications 
for schemes that rely on this interaction in order to enact 
universal quantum computation and/or decoherence 

COIU - ro l 1.2.3.4.5.6.7.8.9 .10.11.19.20.21.22.23.24.25.26.27.28.29.30.31.32.33.34.35.36 .: 

The pure-spin approach cannot be applied in these 
cases, i.e., one must be able to resolve spins in order for 
Heisenberg-based quantum computation to work. 

However, this conclusion does not mean that the 
Heisenberg interaction plays no role in pure-spin dynam- 
ics. In general, this dynamics is non-unitary and as we 
shall see below the constant Jh, characteristic of the 
magnitude of the Heisenberg interaction, appears in both 
the non-unitary term and in the unitary Lam6-shift en- 
ergy, in the p(t) dynamics described by Eq. (|49[) . 

In the general case, in order to obtain a correct effective 
spin Hamiltonian one needs to analyze the exact Eq. (|26p . 
Cts is seen to be a bilinear matrix function of the a s (t) 
and at(t) amplitudes entering as a linear combination 
of their cross-products, while p(t) is a quadratic matrix 
function of the a s (t) and at{i) amplitudes entering in 
separate combinations. Our goal is to express C ts as a 
linear matrix function of p: C t s = As [/>(£)] j so that the 
total initial state information recorded in the a s (t) and 
at(t) amplitudes is compressed into p. To do so we first 
consider the relationship between the amplitudes at the 
initial time t = and t using the time-evolution operator 
Eq. ~ 



a s (t) = U ss {t)a s {0) 
Ot(t) = U ts (t)a 8 (0) 



U st {t)a t {0) 
U tt (t)a t (0) 



(39) 



Evidently, in general it is not possible to establish a 
one-to-one correspondence between the a s (t) and at(t) 
amplitudes. There should be some correlation between 
the amplitudes at t — 0. Let us assume, for example, that 
initial conditions are set up such that we have a linear 
relation between the amplitudes, given by a correlation 
matrix R m (0): 



0,(0) = i? m (o) at (o). 

Then, from Eqs. ((3U) and (gDJl we obtain 
a s (t) = R m {t)a t {t), 

where 



(40) 



(41) 



U H {t) = exp (-iX H (t)si ■ s 2 ) 



(37) 



Rm(t) = (i)^R m (0) + (£)] [c/ f s (t) i? m (0) + 1/" (t)] - 1 . 

(42) 

We will refer to this as the mixed case because the spin- 
density matrix will have both singlet and triplet compo- 
nents if both a s (0) and a t (0) are non-zero. 
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Two other special cases where it is possible to establish 
a one-to-one correspondence between the a s (t) and a t {t) 
amplitudes are: (i) singlet initial state, a t (0) = and (ii) 
triplet initial state, a s (0) = 0. In the singlet (R s ) and 
triplet (R t ) cases, we have respectively 



a s (t) = R Sjt (t)a t (t), 



(43) 



where 



U ss (t)(U ts (t)) singlet case ( , 
U st (t)(U tt (t)) I" 1 triplet case ' 



The triplet case is a particular case of the mixed one 
when R m (0)a t (0) = 0. Assuming the existence of the 
inverse operators (below we consider this issue in detail), 
the matrix F a , a € {s, t, m} in Eq. (|28[) can be rewritten 
as 



F a = Q a (t)\a t (t)) (at(t)\ 



(45) 



where Q a (t) — H fs R a (t). In order to separate unitary 
from non-unitary evolution, we resolve Q a into Hermitian 
and anti-Hermitian parts: 



Pot Qa ~t~ Qa 
Da Qa Qa 



Then, using the identities 



i k k,. k,k; <V,,S. 

K iP K| = (\a t (t)) (a t (t)|V.S, 



(46) 



(47) 



we obtain 



4? = - l -{Pa, P } + lJ2( w *y 



2K lP K\ - P K]K, - K]K lP 



(48) 



where P„ = J2ij (Pa)^ Tij- The first term in Eq, 
describes the effect the singlet-triplet coupling has on the 
unitary part of the system (spin) dynamics, and "renor- 
malizes" the system Hamiltonian (an analog of the Lamb 
shift), while the second one - proportional to the matrix 
D a - is responsible for non-unitary effects in the spin dy- 
namics. Including this Lamb shift into the effective spin 

H" 



Hamiltonian H" 
Eq. (OU 



as 



dp(t) 
dt 



iP Q , we rewrite 



Xij 



([K i)P (t)Kt] + [Kip(t),K5]) (49) 



where Xij = (iD<x)ji- 

Eq. (|4T)|) is one of our central results. It is parametrized 
by the time-dependent functions H*j(t) and Xij(t) which 
in turn can be expressed in terms of matrix elements 



of the Hamiltonian H [Eq. (|22f ]. Since these functions 
depend on which state (singlet/triplet/mixed) is specified 
as the initial state, p{0), we further consider the first two 
cases separately. 

We restrict ourselves to the case of a time-independent 
Hamiltonian. The general case can in principle be re- 
duced to the time-independent one by dividing up the 
time interval into small subintervals and then approxi- 
mating the Hamiltonian by average ones over each time- 
subinterval. 

First we consider the eigenvalue problem 



\e s k) 
\e t k) 



= £k 



\e s k) 
|e t fc) 



Using the closure relation, one obtains 

l (Qe,t - Ql 



T 

Xs,t 



k = !,■■■ ,6. 

(50) 

(51) 



Ps,t — Qs,t + Q S ; 



where 



Qt = {(J2 k ex P(~ i£ kt)H ts \ e sk) (e t fe|) x 

(Efc exp(-is k t) \e tk ) (e^lT 1 } , ^ 

Qs = {(J2k exp(~ie k t)H ts \e sk ) (e sfe |) x 

(Efe exp(-ie fe t) \e tk ) (e sfc |) _1 | , ^ 

are exact expressions in terms of exact solutions of the 
eigenvalue (eigenvector) problem Eq. (|50"|) . In general, 
Eq. ([50]) can be solved numerically using standard Her- 
mitian symmetric eigenvalue problem routines (such as 
in the LAPACK library). Notice the oscillatory behavior 
predicted by Eqs. ((52]) and (j53j) . 

Let us consider an approximate solution to the eigen- 
value problem assuming that H st is a small perturba- 
tion. To zeroth order (H st = 0), the eigenvalue problem 
Eq. ([50)) is reduced to the separate singlet and triplet 
subspace eigenvalue problems 



4V 



H 



and 



( 

,(°)\ - 



eife)=0, = 1,2,3, 

To){ 54 
e° ) ) = ! k = 4, 5, 6 



,(0)\ _ 



|0), k = 1,2,3 



sk/ - |0>, k = 4,5,6, 

where |0) is a zero vector. Up to first order in H st , one 
obtains 





eg) 


4?) = 


(4 0) / 











-1 

Jtt \ ffts 



'sk / i 



,(°)\ 
'tk 



k = 1,2,3, 

k = 4,5,6 
(55) 
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A. Triplet initial state 

After substitution of the vectors Eq. (|55|) into Qt 
Eq. (|52|) we have 



3 6 



Q* = EE 



exp(— ie^t) — exp(— isf't) 



\e£t { 2\H« 



(o) _ _(o) 
k b l 

(°)\ /J°) 



H ts x 



e u \ e ti 



E exp ( _i 4 0) *) 



\l=4 
3 6 



,(°)\ /J°) 



e tl \ e tl 



(56) 



£ £ «p(-fe<°>t)- 

fc=l i,m=4 











(4 0) -4 0) ) 


(4 0) -^) 



(57) 



Neglecting the splitting of triplet energies, which is 
mainly due to the small magnetic Zeeman interaction, 

one can replace 4 = £ t, k — 4,5,6; e t = |X)fe=4 e i°^ 
Also, due to the strong interelectron repulsion in dou- 
ble occupancy states the difference of energies in doubly 
occupied and triplet states is much larger than the differ- 
ence of energies in the singlet, singly occupied state and 
the triplet ones (denote this energy difference, the Heisen- 
berg exchange interaction constant, by Jh — St — e s ) 
so that we can safely neglect the contributions from the 
doubly occupied states. Within these approximations, 
Eq. ([57]) and (|5T|l can be reduced to 



\ 1 — cxpUJat) 

«- Jh 

T _ o sin(J H t) jjts 
Jh 



where 



xi 

Pi 

,(°)\ 



2 1 — CQS( Jfft) jjtS 

Jh 




(58) 



-s2,3 



sl J is the singlet, singly occupied state ( 
are correspondingly doubly occupied states). 

Note that xf 1S n °t necessarily a positive definite ma- 
trix. Indeed, since 



,(()) 



is a positive Hermitian matrix, the sign of the oscillatory 
function g(t) = 2 sin{J ut) / J h determines whether xf 1S 
positive (when g(t) > 0) or negative (g(t) < 0). Addi- 
tionally, this conclusion is supported by the fact that one 
could in principle evolve the state in the backward time 
direction. In this case, g(t) < even at small t. 

Alternatively, and for consistency, one can use the di- 
rect asymptotic expansions at small t in order to derive 
Eqs. (EE|): 



u st {t) = u st (o) + ii st (o)t + t) st {oYj + 

U u {t) = U tt (0) + U u (0)t + J7 t4 (0) V + • 



(59) 



Obtaining from Eq. (|22|) explicit expressions for the time 
derivatives, Eqs. 1|59|) take the form 



U st (t) = - 
U u (t) = I 



iH st t - 
- iH u t 



(H ss H st + 
- (H ts H st 



St TTtt\ t 



H st H 



(JT«) 2 ) f 



(60) 

Inserting these asymptotic expressions into Eq. (|5I|) . we 
obtain 

X f= 2H ts H st t, 
P t = (H ts H st H u - H ts H ss H st ) t 2 . 



(61) 



It can be easily verified that the asymptotic expansions 
(IBTj) coincide with the corresponding expressions (|5"5|) at 
small t (lifting some approximations made above con- 
cerning singlet and triplet energies). 

Also, using Eq. (|6TJ)) one can calculate the triplet states 
population at small t: 

(a t (t)\a t (t)) = (a t (0)|(f/"W) t t/ tt (i)|a t (0)) = l-a t i 2 , 
a t = (at(0)\H ts (H ts ) f \a t (0)) (62) 

Observe that the triplet state probability decreases with 
time, while the singlet state probability, (a s (t) \a„(t)) = 
1 — (a t (t) \at(t)) , increases due to the singlet-triplet cou- 
pling, with the constant at being quadratically depen- 
dent on the singlet-triplet interaction matrix H tB . In 
other words, the total spin of the system is not conserved 
if there is a non-zero coupling between singlet and triplet 
states. 



B. Purity for a triplet initial state 

The time dependence of the purity p(t) = Tr p 2 (t) may 
serve as a measure of non-unitarity in dynamics [p(t) = 1 
iff the state is pure; p(t) is constant iff the dynamics is 
unitary]. With Eq. ([62| we have 

p(t) = Ti-p 2 {t) = 1-2 (ot(t) \a t (t)) + 2 (a t (t) \a t (t)) 2 
= \-2a t t 2 (63) 

at small t. Observe that non-unitary effects are pro- 
portional to the constant 2at which can also be defined 
as the expectation value of the time derivative of the 
matrix xi in Eq. |[6"T]): 2a t = Tr[xf(t)p(0)], where 
Xf = J2ijXuj(t)Tij. In Part II 41 we investigate nu- 
merically how the triplet states population, purity, and 
Lamfr-shift depend on time and other relevant physical 
parameters. 



C. Singlet initial state 

This case differs qualitatively from the triplet one by 
the fact that the inverse operator in Q s Eq. (|5"3"]) is sin- 
gular at t = 0. From the asymptotic expansions 



U ts (t) = -iH ts t - (H ts H ss + 
U ss (t) =1- iH ss t - (H st H ts 



H u H ts ) f 
f (H ss ) 2 ) | 



(64) 
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at small t one finds 

Q s = \H ts (i - iH ss t - (H st H ts + (H ss ) 2 ) 4) x 
(H ts - i {H ts H ss + H u H ts ) iy 1 

(65) 

Assuming that the matrix H ts is non-singular, Q s = jl 
as t — ► 0. However, since Q s is defined on the vectors 
\a t (t)), it operates on the outer product \a t (t)) (a t (t)\ 
in F s Eq. gS]) and \a t (t)) = -iH ts t\a s (0)) at small t. 
Therefore, their product, F s = itH u |a s (Q)) (a s (0)| H st , 
goes to zero as t — > 0. If _ff ts is a singular matrix and 
K(0)) belongs to the kernel of H ts then Q s = 0(t- 2 ) but 
\a t (t)) = - (H ts H ss + H u H ts ) ^ |a s (0)) and F s = 0(t 2 ) 
as t — » 0. When |a s (0)) does not belong to the kernel 
of H ts , the pseudo-inverse of il* s is defined in the Moor- 
Penrose sense (see, e.g., Ref. l49l and references therein). 
Thus, in spite of the singular behavior of Q s at zero we 
have regular behavior of all the corresponding terms in 
Eq. (|49)l and the operator Q s is well defined at t ^ (it 
is not defined in the pathological case where H ts = but 
in this case there is no connection between singlet and 
triplet subspaces and Ct s = 0). 

Keeping up to the second order terms, Eq. (|65|) can be 
rewritten as 

Qs = \l-\ (H u - H ss ^j + | (H ss - H ts H st ^j (66) 

where H ss — H ts H ss (H ts ) 1 is the similarity trans- 
formed matrix H ss and one obtains 

x T s - -f-T + l (h ss -H ss ^ - 

I ^H ss + £P s t - 2H ts H st ^j , 

P s = —H tt + | (H ss + H ss ^ + f (H ss - H ss ^ . 

(67) 

Comparing Eqs. (fBTj) and (|57|) . notice the qualitative 
difference between singlet and triplet cases. While in 
the triplet case the corresponding matrices are basically 
of the second order in the H ts interaction coupling, the 
singlet state matrices, xT an( i ^s> are °f ^he order of 1/t, 
H u , and H ss , respectively, at small t. Observe that at 
small t, is positive when t < and negative when 
t > 0. 

Similarly to the triplet case above, we compute the 
triplet states population 

(a t (t)\a t (t)) = (a. (0)| (U ts (t)) f U ts (t) \a s (0)) =a s t 2 , 

a s = (a s (0)\H st (H st y\a s (0))- (68) 

Here the constant a s is the average of the H st (H st ) 
interaction operator over the initial state |a s (0)). 

D. Purity for a singlet initial state 

Now we have for the purity 

p(t) = Tr p 2 (t) = 1 - 2a s t 2 . (69) 



Alternatively, the rate of non-unitarity, 2a s , can be 
defined as the expectation value of the time derivative of 
the matrix x^ , now in the triplet states \at(t)), since xl 
is defined on \a t (t)) and a t {0) = 0: 

2a s = (a t (t)\ X T S \<H{t)) = Tr \g(t)p(t)], (70) 

where \a t (t)) = -iH ts t\a s (0)) and = {2/t 2 )I, 
£(t)=J2ijXl j (t)T ij at small t. 

IV. THE TIME-CONVOLUTIONLESS FORM OF 
THE MASTER EQUATION 

In the previous section we have examined the possibil- 
ity of constructing a master equation based on the ex- 
plicit structure of the spin density matrix in terms of the 
a St t{t) amplitudes. We found that a meaningful quasi- 
closed form of the master equation can be derived if the 
initial state is assumed to be in a product form of orbital 
and spin functions, correspondingly in singlet and triplet 
states, or in a correlated mixed state. The quasi-closed 
character of the equations we obtained stems from the 
fact that the matrix operators in Eq. (|49|) are defined 
differently in the singlet and triplet cases, although they 
do not depend on the initial state amplitudes, a s (0) and 
O((0) correspondingly. In the correlated mixed case we 
have a correlation matrix R m (0) fixing the relation be- 
tween a s (0) and at(0) amplitudes, and this correlation 
information is present in the operators appearing in the 
dynamical equation [cf. Eqs. (|40f - (l49j) ]. 

In the present section we approach the problem of con- 
structing the dynamical equation using projection oper- 
ator techniques, as applied in the derivation of the TCL 
master equation^ 

A. Derivation of the TCL master equation 

For the analysis of the general mixed case we start with 
the exact equation for p to t 

Ptot = -i [H, ptot] = -iLptot (71) 

where L = L* is the Liouvillean superoperator corre- 
sponding to H. The density matrix p to t contains both 
relevant and irrelevant information about the coordinate 
dependence (bath), which will be averaged out after in- 
tegration over the coordinates: p = Tr^ p to t- Let us as- 
sume that the irrelevant information can be eliminated by 
virtue of a time-independent projection operator P, and 
the relevant information is assumed to be given by the 
projected density matrix Pptot- Using the projection op- 
erator technique ) 18 ! 50 ' 51 one obtains the well-known TCL 
master equation for Pptot 

Pp to tW = K{t)P Ptot {t) + I(t)Q Ptot (0), (72) 

K{t) = -iPL6(t), 

I(t) = K(t)cxp(-iQLt), 

9(t) = (P + expi-iQL^QexpiiLt))^ 1 , 
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where Q = I — P. For simplicity, we restricted ourselves 
to the case where the Liouvillian is time-independent. 
This is an exact, inhomogeneous, first order linear differ- 
ential equation. Both the TCL generator K(t) of the lin- 
ear part and the inhomogeneity I(t) are explicitly time- 
dependent superoperators which do not depend on which 
initial state ptot(0) is taken. Observe that the inhomoge- 
neous term vanishes if the initial state satisfies the rela- 
tion Qptot(O) = 0, i.e., if 



Pptot(O) = ptot(0). 



(73) 



This coincides with having factorized initial conditions. 18 
Below, we consider when this is possible in our setting. 

The condition that Pptot contains all the relevant in- 
formation means that the following condition 



p = Tr^ptot = Tr^Ppt, 



(74) 



together with P 2 — P should be imposed on P as a 
projection operator. Let us take any orthogonal decom- 
position of the unit operator Ip on the coordinate state 
space span {f ai , f t } i=1 2 3 



If = J2\fsi) (U + lft) (ft\, 



(75) 



i.e., a collection of projection operators LT a that satisfy 

n a n fc = s ab n b , ^ n a = i P . (76) 



Then, following Ref. [l8| we can define a linear map by 
means of 



(77) 



where N a — Trf {n a }. It is easy to check that this su- 
peroperator fulfills the above conditions. As an example, 
let us consider the case where IT a = n s = ^2 i=1 \ f S i) (fsi\ 
and — LT t = \f t ) (ft\. Then, the six-by-six density 
matrix 



p ss p st 

/>!■■' — 1 pts ptt 



(78) 



is transformed under the map Eq. ([77)) to 

Pptot = ( P'f ) (79) 

where — \Y^i=iPii- Notice, however, that the 
projection superoperator P must be realizable as a su- 
peroperator independent of ptot, i-e., the matrix state 
Eq. (|79[) must be obtained by multiplying the matrix 
(175)) from the left and right by some matrices inde- 
pendent of ptot- It is clear that by taking Tl a = 



{U si = \ f si ) (f si \ , i = 1, 2, 3; II t } the projection superop- 
erator P can be written as 



Pptot = 



diag(p-) 







r 



— ^ E si p^t E S i + E t ptotE t 
i=i 

(80) 

where diag (pff ) is a diagonal matrix with pff, i — 1, 2, 3 
entries on the diagonal; E S i and E t are orthogonal pro- 
jection matrices independent of p to t and defined by 



Esi = 



U 




E t = 




/ 



(81) 



Here, li is a diagonal matrix with all zero entries except 
1 in the i-th entry. However, without changing matrix di- 
mension it is not possible to change the diagonal entries 
p\l in Ea.pO]) to p™ by multiplying the matrix Eq. (j80|) 
by some diagonal matrices independent of ptot- Hence, 
in order to comply with the requirement of state inde- 
pendency one should extend the definition of the projec- 
tion superoperator to Eq. (|80f . Observe that the state 
Eq. ([80)) contains two extra pieces of information: be- 
sides the probability of single state occupation, pf f , also, 
the probabilities of double state occupancy, p 2 2 an d /°If> 
which will be compressed after averaging Tr^Pp t ot to 
p£. Then, it follows from Eq. (80)) that 



Qptot — E s ptotE t + E t ptotE s 



&3 



E si Ptot E, 



•V 



(82) 



where E s = J2i Esi- 

In principle, one can "compress" diag (pf, 
by applying a rectangular 6-by-4 matrix 



C 




in Eq. Mf\ 



(83) 



so that 
Pptot 



C T EsiPtotEgi + E t ptotE t j C 



E t pf! o 



1x3 



0, 



3x1 



(84) 



However, the transformation Eq. (|84|) is not a projec- 
tion operation: C T C C. 

In numerical computations it is often convenient to 
rewrite matrix operator equations in tensor product form 
(see for example Ref. and references therein) 

Y = CXB T <= 



\Y) = (B®C)\X) 



(85) 



where Y, C, X, B e C" xn , B <g> C e C" x ™ are matrices 
and we introduced the Liouville vector space \X) G C n 
by applying the "vectorization" operation to a matrix: 



X 



\X) 



x. n 



(86) 
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which amounts to a "stacking" of matrix columns. It is 
easy to check that (X \Y) = Tr (X^Y). 

Using these rules, Eqs. ([85]) and ([86]) . we can rewrite 
the superoperators in Eq. (|72p in the tensor product form 

3 

L = I®H-H T ®I, P = Y^E s i® E si + E t ® E t , 
Q = E s ®E t + E t ®E s + J2E si ® E sj . 

(87) 

Then, with the use of identities (A gfl)(C®D) = AC ® 
BD and (A + B)®C = A®C + B®C, one obtains 

QLQ = ^2 (E t ® E si HE sj - E si B T E sl ® E t ) + 

ij 

(E s ® E t HE t - E t B T E t ® E.) + 
i E « ® E sjHE sk - 

E S jH T E s k ® + 

(E t H T E SJ + E S3 H T E t ) ®E sl }, (88) 



PLP = J! ® - E sl H T E sl ® E sl ) + 

i 

E t ® E t HE t -E t H T E t ®E tl (89) 

QLP = E t ® E s HE t - E s H T E t ® E t + 

Y {Esi ® E t HE si - E t H T E sl ® E si ) + 

i 

^2 (Esi ® E S jHE S i — E S iH T E s j ® E S j) , 

PLQ = (QLP)t. (90) 

Having obtained the basic superoperators in tensor 
product form, in the next step we can consider how to 
compute the superoperators Kit) and I(t) in Eq. (j?2"|) 
and then find the conditions under which Eq. (|73p 
(factorized initial conditions) is fulfilled. Inspection of 
Eq. I|72p shows that in order to calculate K(t) and I(t) 
one needs to calculate the exponentials exp(— iQLQt) 
and exp(iLi) and take an inverse of a matrix to obtain 
9{t). Calculation of exp(iLi) does not introduce any dif- 
ficulties and can be reduced to the eigenvalue problem 
Eq. (|50[) . namely, the eigenvalue problem for the Liouvil- 
lean superoperator 

(I®H-H T ®I) \E nm ) = e nm \E nm ), (91) 
n,m = 1, 2, • • • , 6 



is solved in terms of eigenvalue solutions to Eq. (I50[) : 

£ nm = £n - £ m and E mn = \e n ) (e m \ where 




and \E nm ) is obtained according to the rule Eq. 
Then, 

exp(iLt) = Y exp(ze„ m i) \E nm ) (E nm \ . (92) 

nm 

A direct calculation of exp(— iQLQt) does not seem to 
have a similar simple solution due to the rather compli- 
cated tensor structure of QLQ Eq. ((88)) . However, one 
can overcome this problem by simplifying 9(t) to 

6(t) = (cxp(-iQLQt) [exp(iQLQt)P + Q cxp(iLi)]) -1 
= (P + Qexp(iLt)) -1 exp(iQLQt), 

(93) 

where we have used the identity exp(iQLQt)P — P and 
unitarity of exp(— iQLQt). Since the 6{t) superoperator 
in K{t) operates on the state P |ptot) and the Q-projected 
exponential is cancelled in the I(t) superoperator, the 
TCL equation can be rewritten in the form 

J t \ p P^{t)) = -iK{t){\P Ptot {t)) + \Q Ptot {Q)))m) 
K{t) = PL(P + Qexj>(iLt)y 1 (95) 

free from exp(— iQLQt). In fact, as one can see K(t) and 
I(t) are the same superoperators. 

Eq. along with (prS)) are our central results. Know- 
ing how to calculate the exponential Eq. (f9"2"|) . it is then 
straightforward to compute Kit) Eq. (|95|) by applying 
standard matrix multiplication and inversion routines. 



B. Separation of unitary and non-unitary dynamics 

Condition (|73| for a factorized initial state can be 
rewritten as 

Q^0)=(*" m -;^ m ^)=0. (96) 

It can be fulfilled exactly in the following two cases: (i) 
the singlet case: Ot(0) = and any of the singlet state 
amplitudes a s j(0) ^ with all the others being zero, 
a S j(0) = 0, j ^ i; (ii) the triplet case: a s (0) = 0. In 
other cases, including the mixed one, Eq. (|94[) will not 
be in closed form and will contain a non-zero inhomoge- 
neous term. Notice that unlike the Lindblad formulation 
[Eq. HIS}], the superoperator K(t) Eq. |S5} does not de- 
pend on which initial state - singlet or triplet - is taken. 

Henceforth we neglect the inhomogeneity for simplic- 
ity. To separate unitary effects from non-unitary ones, we 
decompose the dynamics generator K{t) into Hermitian 
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and non-Hermitian parts: 



D. Purity 



K(t) = K+{t)+K_{t), 
K ± {t) = \ {K{t)±K\t)} 



(97) 



Transforming to the interaction representation (denoted 
by a superscript hat) we have for the unitarily trans- 
formed state 

|PptotW) = exp (i J K+(r) dr^j \Pp to t(t)) (98) 

the following equation of motion: 
d 



dt 



\P(>tot(t)) = -iK_{t) \P P tot(t)) , (99) 



where 



K-(t) = T^cxp U / K + (r)dT 



K-{t)T^exp {-i J q K+(T)dT^j . (100) 

Here, (TU)T<_ denotes the (anti)-chronological ordering 
operator. One can see that Eq. (|99p describes non- 
unitary dynamics since the generator K^_(t) = —K~{t) 
is an anti-Hermitian superoperator. If K—(t) = 0, then 
the evolution is unitary 

|PAot(*)> =r<_exp (-ij K + (T)dr^j |Pp tot (0)). 

(101) 



C. Short time expansion 

Expanding K(t) in powers of t, one obtains 
K+{t) = PLP + ({{QLP)\QLP),PLP} - 

2{QLP)\QLQ){QLP)) t - + --- , 
K-(t) = -i{QLP)\QLP)t + 

[{QLP)\QLP) 1 PLP] t - + --- (102) 

where {•,•} denotes the anti-commutator. In the zeroth 
order approximation K(t) = PLP, Eq. (|94[) is reduced 
to a unitary matrix equation of the form 



/<;; * = i,2,3 

p u = -i[H tt ,p tt ] 



(103) 



which describes a unitary evolution of the triplet state 
density matrix p a under the Hamiltonian H tt with the 
singlet state population probabilities pff (t) = pff (0) be- 
ing constant in time. Let us consider the non-unitary 
effects which are induced by the first term proportional 
to t in the expansion of K _ (t) Eq. (|102[) . 



Observe that Tr Pp tot = Tip but Tr(Pp tot ) 2 ^ Tip 2 . 
In order to obtain p from Pptot we first have to apply 
the "compression" transformation Eq. (|8"4")) so that the 
purity becomes 

pit) = Tvp 2 {t) = (P Ptot (t)\ CC T P Ptot (t)CC T ) (104) 

where 



CC 



T _ ( -hx3 
/3X3 



J; 



3x3 




(105) 



Then, using Eq. (|94[) we obtain an expression for the time 
derivative 



dt 



Trp 2 (t) = i(P Ptot {t)\{K\t)CC J 



CC T K(t))P Ptot (t)CC T ) . (106) 



It is easy to check that the first term, PLP, in the ex- 
pansion of Kit) does not contribute to Eq. (| 1 06|) and the 
first non-zero contribution comes from the first term in 
the expansion of K_(t) [Eq. (|102p ]. which can be written 

as 



dt 



Trp 2 (t) = -2t (QLPp to t(0)| QLPCC T Pp tot (0)CC T ) . 

(107) 



Further, with the use of Eq. (|90p one obtains 



QLPp tot (0) 



Cll C12 

-4 



(108) 



Cll = [P ss ,diag(^(0))], 

ci2 = H st p tt (Q)-dia,g(p ss (0))H st , 



and 

QLPCC T P Ptot {0)CC T 



(109) 



di 2 
dai 

di2 = H st p u (0)-p s a s (0)H st , 
d 21 = H ts p s a s (0)-p tt (0)H ts , 

where p s a s (0) = EiPf/(0)^3x3 and Eq. CTO7|) is reduced 
to 

jTrp 2 (t) = -2t{Tr(P*Vr(0)-p"(0)P ts ) x 

(diag (p ss (0))H st - H 3t p tt (Q)) + 
Tr(iPV*(0) -p s a s (0)H st ) x 
(p tt (0)H ts ~ if ts diag(p ss (0)))} (110) 

Finally, from Eq. (|110[) one obtains 
( (H st H ts ) u , 

±Trp 2 (t) = -At f^t C rrt%$ {0 li ° ; (HI) 
dt F K ' I (a t (0)\H ts H st |a t (0)) , v ; 

[ triplet case: p"(0) 7^ 0, 
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which is equivalent to the results we obtained in the Lind- 
blad analysis for the purity, Eqs. ([55)) and ([53]) respec- 
tively. We provide a detailed discussion of the behavior 
of the purity and other physically relevant quantities in 
Part 1141 



E. TCL vs Lindblad-type dynamics 

Note that the Lindblad-type and TCL Eqs. ([49]) and 
(|94[) are not mathematically equivalent formulations. In 
the Lindblad-type formulation we allowed correlations 
between initial singlet and triplet states to be incorpo- 
rated into the dynamics, whereas in the TCL formulation 
we assumed that all operators appearing in the dynam- 
ical equations are totally state independent. The first 
approach gave us the flexibility needed to cover mixed 
states, whereas the TCL equations become non-closed 
[inhomogeneous term in (|94|) is non-zero] in the mixed 
state case. 



V. SUMMARY, DISCUSSION, AND 
CONCLUSIONS 

In this work we have shown that the time-evolution of 
the spin-density matrix p(t) that describes a system of 
two electron spins is in general non-unitary when spatial 
degrees of freedom are accounted for, without coupling 
to a "true" bath. The non-unitary effects (e.g., the to- 
tal spin is unconserved) in the pure-spin evolution are 
due to a non-zero coupling between singlet and triplet 
states. Our primary focus in this work was the deriva- 
tion of dynamical equations for the spin-density matrix, 
and a corrsponding study of the non-unitary effects. In- 
voking standard ideas from the theory of open quantum 
systems, one can define a "system + bath" by formally 
associating the spin variables with the "system" degrees 
of freedom, and the spatial variables with the "bath" de- 
grees of freedom. With system and bath thus defined, 
one can apply the standard machinery of open quantum 
systems theory to attack the above problem. We did 
this by first constructing an analog of the Kraus opera- 
tor sum representation for p(t) and then deriving master 
equations for p(t) in the Lindblad and TCL forms. 

The Kraus representation [Eq. (fl5[) ] has the shortcom- 
ing that its Kraus operators are dependent on the initial 
state p(0). This is a consequence of the fact that an arbi- 
trary initial state p(0) cannot be represented as a product 
state with some reference state in coordinate space due 
to the Pauli antisymmetry principle, which is imposed on 
the total wavefunction. As a result of this dependence, 
the complete positivity property of the mapping (|15p is 
not guaranteed. 

The Lindblad-type master equation Eq. ([49]) describes 
the exact dynamics of p(t) (it is not a Markovian ap- 
proximation in our case). Again, the generators of this 
dynamics are not totally independent of initial condi- 



tions. This is to be expected; after all, the p-dynamics 
is inherited from the total system unitary dynamics de- 
scribed by Eq. (jTTj) . In the total Hilbert space, the 
state is defined by 11 real parameters [6 complex ampli- 
tudes {a s (t), a t (t)} minus the normalization condition], 
while p is defined by 5 real parameters ( "system degrees 
of freedom"): 3 amplitude moduli |a«(i)|, i = 1,2,3, 
and two relative phase angles between the amplitudes 
a t i(t), a t2 (t) and a n (t), a t3 {t) \ftuj{t) = a tl {t)a* ti {t) , 
p s (t) = 1 — J2i Ptii(t)]- Thus, we have 6 extra degrees of 
freedom [4 real parameters defining two complex double- 
ocupancy amplitudes, a S 2(t) and a s ^(t), and two phases 
of the amplitudes a s i(t) and ati(t)] in the reduction 
^tot - > P, which should be present in the dynamics equa- 
tion if it is exact. The three complex-valued equations 
(|4T))) . that establish relationships between amplitudes a s 
and at at the initial moment, allow us to express 6 ex- 
tra degrees of freedom as a function of system degrees of 
freedom. Fixing the correlation matrix R m (0) appearing 
there, defines a domain in the space of 11 real param- 
eters so that dynamical operators in Eq. (|49p acquire a 
dependence on the domain on which a correlated mixed 
state is defined. Of course, during time evolution R m (t) 
evolves according to Eq. (|42[) and the {a s (t), a t (t)} point 
can move out of the domain, fixed by R m (0), on which 
the initial states were defined. In this case, R m {t) will 
define a new domain at time t, which can be considered 
as the initial state domain for later times. 

It turns out that the Heisenberg interaction does not 
affect the pure-spin dynamics as long as the coupling be- 
tween singlet and triplet states is neglected - see Eq. ([3"5[) . 
If there is a non-zero coupling, the Heisenberg interac- 
tion modifies the corresponding non-unitary terms: the 
Heisenberg interaction constant Jh appears in Eq. (|58p . 
It also appears in the unitary Lamb shift term [the ma- 
trix P t , which defines this shift, depends on Jh]- Note 
that the matrix functions in Eq. ([55[) . responsible for 
the singlet-triplet coupling, depend quadratically on the 
magnitude of the interaction between singlet and triplet 
states, described by the H ts matrix elements. 

The TCL Eq. ([9~4")> is also exact and describes a non- 
unitary evolution; the dynamical generator K{t) is non- 
Hermitian and does not depend on initial conditions. 
However, the TCL equation is not closed with respect 
to the \Pptot(t)) state since it contains an inhomoge- 
neous term, which is non-zero in the mixed initial state 
case. In the non-mixed case defined by Eq. (j96[) . when 
the TCL equation becomes closed, we found a short- 
time expansion [Eq. (|111[) ] for the purity function p(t), 
whose time-dependence signifies non-unitary effects in 
time-evolution. For consistency, we obtained the same 
short-time formulas within the Lindblad-type formula- 
tion, Eqs. (|63p and (|69p . In Part II,— we provide 
a detailed numerical example demonstrating these non- 
unitary effects. 

In conclusion, the spin-density matrix completely de- 
scribes the pure-spin dynamics of a two-electron system 
in the case where the electrons cannot be spatially re- 
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solved. The formalism we have developed can be gen- 
eralized in several directions. First, one could take into 
consideration excited-state orbitals, which will increase 
the dimensionality of the total Hilbert space. Another 
interesting generalization is to consider several interact- 
ing QDs along the lines developed in Refs. Isllollloll and 
to derive for this few-body case the corresponding spin- 
density matrix dynamics. 
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